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Abstract 

We present a systematic study of different methods for analytic calculations of diffuse 
spectra of ultra-high energy nuclei. Nuclei propagating in the intergalactic space are 
photo-disintegrated and decrease their Lorentz factor due to the interaction with Cosmic 
Microwave Background (CMB) and Extragalactic Background Light (EBL). We calculate 
the evolution trajectories in the backward time, that describe how atomic mass number 
A and the Lorentz factor V change with redshift z. Three methods of spectra calculations 
are investigated and compared: (i) trajectory method, (ii) kinetic equation combined 
with trajectory calculations and (Hi) coupled kinetic equations. We believe that these 
three methods exhaust at least the principal possibilities for any analytic solution of the 
problem. In the most straightforward method (i) only trajectory calculations are used to 
connect the observed nuclei flux with the production rate of primary (accelerated) nuclei 
Aq. In the second method (ii) the flux (space density) of primary nuclei, and secondary 
nuclei and protons are calculated with the help of kinetic equation and trajectories are 
used only to determine the generation rates of these nuclei. The third method (Hi) consists 
in solving the complete set of coupled kinetic equations, written starting with primary 
nuclei A$, then for Aq — \ etc down to the A of interest. The solution of the preceding 
equation gives the generation rate for the one which follows. An important element of the 
calculations for all methods is the systematic use of Lorentz factor instead of energy. We 
consider here the interaction of nuclei only with CMB, this case is particularly suitable 
for understanding of the physical results. In paper (II) of this series the EBL will be also 
included. Estimating the uncertainties of all methods discussed above, we conclude that 
the method of coupled kinetic equations gives the most reliable results. 

1 Introduction 

The observation of particles with ultra high energies (E > 10 18 eV) has a fundamental 
importance in high energy astrophysics and, maybe, more generally in physics. These 
particles, hereafter referred to as Ultra High Energy Cosmic Rays (UHECR), attract 



much attention in the recent years as the most energetic particles ever observed and 
because they are messengers from cosmic accelerators to the highest energies. 

The observational data on UHECR can be divided into three categories: spectra, 
chemical composition and correlation with astrophysical sources. 

Several theoretically well established features are predicted for the energy spectra 
of UHE protons interacting with the Cosmic Microwave Background (CMB) radiation. 
Most noticeably they are: (i) the Greisen-Zatsepin-Kuzmin (GZK) feature [I], a sharp 
steepening of the spectrum at E ~ 5 x 10 19 eV due to the photo-pion production [I]; (ii) 
a rather faint feature, the dip, at energies 1 x 10 18 - 4 x 10 19 eV, caused by the e+e" 
pair production [2], and (hi) an even more faint and sharp peak at energy 6.3 x 10 19 eV 
produced by an interference effect in the proton interaction with CMB |3j. 

Two of these features, the GZK cutoff and the dip, most probably are already de- 
tected. The observation of the GZK feature in the differential spectrum of HiRes [1] is 
strengthened by the measured value of E-^ii in the integral spectrum, which precisely 
coincides with the theoretical prediction. The Auger data [5j EJ El [8] also show a sharp 
steepening of the spectrum, which however does not fit well the theoretical shape of the 
GZK feature, especially in the data of 2010 |8j. 

The dip is very well confirmed by the data of Hires, AGASA and Yakutsk detectors 
[TOl [IT], and by the Auger data of 2007 [12] though with a larger \ 2 - The agreement 
of Auger data 2010 with the dip becomes even worse. 

One may notice that both features discussed above are signatures of a pure proton 
composition. Therefore, agreement of the data with the dip and GZK cutoff, most notice- 
ably in HiRes data, can be considered as an indirect evidence for the proton-dominated 
composition. 

The direct observation of UHECR mass composition is a very difficult task (for a 
review see [13] and [H]). At present, observations are contradictory. While HiRes |15] . 
HiRes-MIA [16] and Yakutsk [T7] detectors favor a proton-dominated flux at energies 
E > 10 18 eV, Fly's Eye [18] , Haverah Park [19J and AGASA [20] and, recently, Auger [UJ 
indicate a mixed nuclei composition . 

The most dramatic conflict exists now between Auger \21\ [22] and HiRes \23\ [24] data. 
Both detectors measure two quantities which characterize the mass composition: A max 
as function of energy and RMS(A max ), the width of distribution over A max . While data 
of HiRes evidence for pure proton composition, the data of Auger show nuclei mass com- 
position, which becomes steadily heavier with increasing energy and approaching almost 
pure iron composition at energy ~ 3.5 x 10 19 eV. On the other hand, the observation of 
the correlation with AGN [25] seem to be an evidence of proton-dominated composition 
starting from energy (5 — 6) x 10 19 eV |26j . 

From a theoretical point of view, there are two models with different predictions for 
the mass composition. The dip model predicts (almost) pure proton composition starting 
from energies E > 10 18 eV [27] : the second model [28] predicts a mixed composition. 
As shown in [6], both models in their present form contradict the Auger data on mass 
composition. However, the mixed model has more power in the explanation of the Auger 
mass composition: due to presence of several free parameters, the predictions can be, in 
principle, adjusted to explain a wide range of different mass compositions. 

The discussion above demonstrates the importance of a theoretical study of nuclei as 
carriers of the UHE signal. The theoretical interest to UHE nuclei is also supported by 
the difficulties in reaching the highest energies by the acceleration of protons in astro- 
physical models. In the case of nuclei this problem is ameliorated because the maximum 



acceleration energy increases by a factor Z, the charge number of a nucleus, which is 26 
for Iron. 

Historically, the interest to UHE nuclei started from the hope to solve the UHECR 
puzzle, i.e. the absence of the GZK cutoff in observations, with the help of nuclei as 
signal carrier. In this case the energy of a CMB photon in the rest system of a nucleus 
with atomic mass number A is A times lower than for a proton of the same energy, thus 
photo-pion production is suppressed. 

The first calculations of the UHE nuclei photo-disintegration by the CMB radiation 
has been performed by Stecker [29] and by Puget and Stecker [30j, where infrared radiation 
was also included. The first calculations of the diffuse nuclei spectrum have been made 
by Berezinsky, Grigorieva, Zatsepin [31] and Hillas [32], both works have been presented 
at the same conference and published in the same volume of the proceedings. In |31] it 
has been demonstrated that the steepening in the heavy-nuclei spectra such as Fe, Al, C 
etc occurs at the Lorentz-factor determined by equality of adiabatic and pair-production 
energy losses, i.e. at energies lower than the GZK cutoff for protons. 

Propagating in intergalactic space, UHE nuclei interact with the CMB and Extragalac- 
tic Background Light (EBL), i.e. Infrared, Visible and Ultra Violet radiation, experiencing 
two main processes: (i) photo-disintegration and (ii) e + e~ pair production, that reduce 
the nucleus kinetic energy. The systematic study of the photo-disintegration started from 
the pioneering works by Stecker [291 GSB E3J EI] • These works are based on a very con- 
venient parameterization of the photo- nuclear cross-sections, widely used in the study of 
UHE nuclei propagation. Recently, a refined parameterization of this cross-section was 
presented in [35J. 

Using the parameterization of the photo-nuclear cross section by Stecker [33], updated 
in 1999 [34J, or using other approximations as in [35], many papers have been published 
in the last decade, all (except [36l|37]) based on a Monte Carlo (MC) approach to the 
study of UHE nuclei propagation in astrophysical backgrounds [281 EH ESI ESI ESI HQl El 

ma us mi hu hu sa mi . 

We present and compare here three new methods for the analytic calculation of UHE 
nuclei spectra: the evolution-trajectory method, kinetic equation method combined with 
trajectory calculations and the method of coupled kinetic equations. In this paper (paper 
I) we include only CMB as the radiation with which UHE nuclei interact. The main 
emphasis of this paper is given to the theoretical issues of these methods. It is more 
convenient to do it with CMB only, because the Planck spectrum and the exact knowledge 
of the cosmological evolution of this radiation make formulae simpler and the description 
of the method more transparent. The phenomenological predictions obtained with CMB 
only are always relevant for the highest energy part of the spectra. 

In a forthcoming paper (paper II) we will include also the interaction with EBL radi- 
ations, and focus more on the observational applications of our new computation scheme. 

The paper is organized as follows: in section [2] we calculate the energy losses of nuclei 
on the CMB radiation, and compute the evolution trajectories. In section E] the combined 
method is presented, where space particle density is obtained from kinetic equation, and 
trajectories are used only for calculations of the generation rates. In section [5] the method 
of coupled kinetic equations is developed. Conclusions are presented in section [5j In 
the appendixes we give a detailed calculation and discussion of the following technical 
problems: generation rates calculated from number of particles conservation, in appendix 
[Al the ratio of energy intervals at different redshifts, in appendix [Bj the analytic solutions 
of kinetic equations, in appendix [Cl and the comparison of the secondary nuclei and 



proton fluxes, in appendix [Dl 



2 Energy losses and evolution trajectories 

Propagating through background radiations, mainly CMB and IR, nuclei decrease their 
Lorentz factor T, due to e + e~ pair production, and atomic number A, due to photo- 
disintegration. Lorentz factor decreases also due to the expansion of the universe (adia- 
batic energy losses). In this section we study the energy losses and the evolution of nuclei 
due to their propagation, considering the evolution trajectories in time (or redshift z) 
along which A and V are changing. We will consider these trajectories in the backward 
time, starting from the point of observation (T, A, zq = 0) and increasing z, so that T and 
A increase too. 

An important ingredient of our method consists in the use of the Lorentz factor of 
the particles instead of their energy. This approach gives many simplifications in the 
theoretical study, in particular the approximate equality of Lorentz factors of all three 
particles participating in photo-disintegration process A — > (A — 1) + N, where N is a 
nucleon. 



2.1 Energy losses on CMB 

In this paper we consider the interaction of nuclei only with CMB. Adiabatic energy losses 
and pair production change only the Lorentz factor and photo-disintegration changes only 
the atomic mass number A. Then for E = TAuin we have 

l_^E_l_dr l_dA 

E~dt ~ f~dt + A~dt' ^ ' 

For any process the energy losses in continuous approximation for nuclei or protons 
interacting with CMB at z = can be written as [50] : 



1 dv t r°° r / f \ 



(2) 



where e is the photon energy in the nucleus or proton rest frame, eo is the threshold of the 
considered reaction, /(e) is the mean fraction of energy in the laboratory system lost by a 
nucleus in a single interaction, i.e. inelasticity, <r(e) is the cross-section and T = 2.726°K 
is the CMB temperature; units h = c = k = 1 are used. 

At redshift z the number of CMB photons is (1 + z) 3 larger than at z = and their 
energies are (1 + z) times higher. Then the energy loss at arbitrary z is given by 

P(T,z) = (l + zfp [(l + z)T} (3) 

In the following we will specify the three processes relevant for our calculations: expansion 
of the universe (adiabatic energy losses), pair production and photo-disintegration. 

(i) Adiabatic energy losses. 
Adiabatic energy losses are given by 

Mz) = -~ = H(z), (4) 

where H(z) = Hqv(1 + z) A Q. m + Q\ is the Hubble parameter at redshift z with Hq = 
72 km/sMpc, 9, m = 0.238 and tt A = 0.716 according to WMAP data [EE]- 



(ii) Electron-positron pair production 

This process can occur if the energy of the background photon is larger than 1 MeV in 
the rest system of the UHE particle. The inelasticity and cross-section for a nucleus are 
simply related to the corresponding quantities for the proton pair-production: 
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where Z is the nucleus charge number, is the nucleon mass, o" pa i r is the pair production 
cross section for protons [i = p) or nuclei (i = A). Using these equations we can rewrite 
Eq. © as 



/^airtr) 



1 alT 
f~dt 



z 2 



A 2vr 2 r 2 



de ^air 



In 1 



exp[ 



2TT 



(6) 



The expression above is valid for z = 0. For arbitrary z one should use Eq. ([3]). In our 
calculations we use the function /3p air (T) for protons as computed in [10]. The quantity 
rr = /3p~ a ^ r has an important physical meaning being the characteristic time of the Lorentz 
factor decreasing. 



(ii) Photo- disintegration. 

The photo-disintegration energy losses on CMB at z = are given by 
_ A fT , s IdA T 1 f°° , . . . 
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where z/(e) is the average multiplicity of the ejected nucleons. 

The quantity /3^ s determines the time-scale of the total photo-disintegration of a nu- 
cleus r tot ~ l//?di s , while T4 = ((M/cii) -1 gives the mean time for one nucleon loss. 

At arbitrary z the pair-production energy losses /3 pa i r (r,z) and photo-disintegration 
energy losses /3^ s (r, z) are given by Eq. (J3]). 

In the present paper we use the remarkable collection of nuclear cross-sections and 
their parametrization from the works by Stecker et al [331 [34] . 

Depending on the photon energy e in the nucleus rest frame, one can distinguish 
two different regimes of photo-disintegration, namely the low-energy regime e < 30 MeV 
and the high-energy regime: 30 < e < 150 MeV. At energies e > 150 MeV the photo- 
disintegration process is not important [33j[3l] and we include at these energies the photo- 
pion production. 

In the low-energy regime the leading process of photo-disintegration is represented by 
one and two nucleon emission, A + 7— s>(A — 1) + A^,A + 7^(A — 2) + 2N , and the 
photo-disintegration cross-section is dominated by the Giant Dipole Resonance (GDR) 
[33], [34] . At higher energies a multi-nucleon emission regime takes place (see [33] [34] ) . 
In this regime the photo-disintegration cross-section can be approximated as constant 
[331 [34] . The main contribution to photo-disintegration is given by the GDR cross-section 
with mean number of emitted nucleons 2 > v > 1. 

We have calculated the energy losses for all nuclei with existing data, presenting some 
of them in Figs. [1] and [2] as function of the Lorentz factor. 

As discussed below, an important quantity that characterizes nuclei propagation and 
energy spectra is the critical Lorentz factor r c , defined from the condition of equality of 



c 



photo-disintegration lifetime ta = (dA/dt) 1 and the characteristic time rp for changing 
of the Lorentz factor, 

ta(T c ) = r r (r c ), (8) 

where tt is given by the sum of adiabatic and pair-production losses as rf 1 = t~ x +T~^ ir . 

Eq. ([8]) is written for z = 0. The same condition for an epoch z determines T c at 
redshift z, for which we will use the notation T c (z). 
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4 x 10 8 


1.2 x 10 9 


4 x 10 9 


3.8 x 10 s 


2.2 x 10 9 


2.1 x 10 9 


2.3 x 10 9 


1.9 x 10 9 



Table 1: Values of r c for some selected nuclei 



The values of T c calculated for some nuclei from Eq. [8] are listed in Table [TJ From 
this Table and Fig. [2] one may notice that r c is almost the same for all heavy nuclei 
with a value around (2 — 2.5) x 10 9 . The values of T c calculated from the condition 
Pa(T c ) = /3r(r c ) (see Figs. Q] and ED differ but little from the values of Tabled) 



2.2 Evolution trajectories 

A primary nucleus Aq accelerated to large Lorentz factor T g at redshift z g is soon photo- 
disintegrated and then evolves as a secondary nucleus with decreasing A and T as redshift 
z decreases. 

Our study of evolution goes in the backward time, whose role is played by redshift. 
We consider as initial state a secondary nucleus A with Lorentz factor T and evolve it 
to larger z with increasing of A(z) and T(z), until A{z) reaches Aq at z g (index g here 
and henceforth implies generation, i.e. acceleration, of nucleus Aq). Evolution of T from 
zq = to z g gives T g of nucleus Aq. This backward-time evolution is governed by two 
coupled differential equations, which for the variable z reads 



l_dA 

A dz 



dt 

dz 



P dis (T,A,z) 



1 dT 

r 



dt 

dz 



/3 pair (r,A,z) + 



l + z' 



being the relation between time and redshift dt/dz 



dt 

dz 



1 



(i+z)H( z y 



H{z) = H Q J{l + zY£l m + £l 1 



(9) 



(10) 



The term 1/(1 + z) in Eq. ([9]) corresponds to adiabatic energy losses dT/dt = —TH(z). 

The numerical solution of the coupled equations Q with initial values A and T at zq 
gives the evolution trajectories A{z) = A(A,T, zq, z) and T(z) = Q(A, T, zq, z), where the 
first three arguments describe the initial condition, in most cases with z$ = 0. 

In principle the solution A(z) is given for continuous A. We shall refer to the cor- 
responding trajectories as A-continuous. In this solution we use /3^ s smoothly changing 
between two integer values of A. 

In the kinetic-equation approach, which is used for calculation of particle density 
ua (I\ z), we assume that A is constant until it reaches due to A-evolution A + 1 (or 
A + 2). The trajectories with this evolution of A we call A-jump trajectories. They are 
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Figure 1: Energy losses for light nuclei due to photo-disintegration and pair production on 
CMB (red full line and blue dotted line, respectively) and adiabatic energy losses given by H 
(green dashed line). 





Figure 3: A sketch of T(z) and A(z) backward-time evolution with T(z) evolution shown by 
continuous curve and integer values of A shown by crosses. A(z) is calculated as continuous 
quantity only for determination of z where integer values of A(z) are reached. The transition 
e.g. from A to A + 1 nuclei (or A + 2), is assumed to occur instantaneously and marked by 
crosses. A nucleus between two crosses is considered as one with fixed A. Evolution of T 
between two crosses is calculated for fixed A. 
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Figure 4: Calculated evolution of the atomic mass number A(z) (left panel) and Lorentz factor 
r(z) (right panel) for a deuterium nucleus observed at z = with various energies as labeled. 



calculated assuming in Eq. Q for T-evolution A = const, until A reaches A + l (or A+2). 
Accordingly, we use in this case the jump behaviour of /3^ s . 
In fact, both methods give practically identical results. 

A calculated trajectory is sketched in Fig. [3]as T(z) = Q(T, A, z) with integer A marked 
by crosses. In the kinetic equation approach we will interpret this trajectory assuming 
that the nucleus A is produced in an instantaneous decay of the nucleus A+l, then it lives 
as nucleus A and finally instantaneously decays to A — 1. The intervals between crosses 
are determined by continuously changing A(z) setting AA = 1 (or AA = 2). This is the 
most natural interpretation, because a nucleus with continuously changing A just does 
not exist in nature, while instantaneous production of nucleus A from decay of A + 1 is 



a realistic photo-disintegration assumption, as well as an assumption that after decay of 
A + 1 in the forward-running time a nucleus A propagates with fixed A until instantaneous 
A — > (A — 1) decay (we always imply below the case of A + 2, too). 

In Fig. d] the calculated evolution trajectories A(A, T, zq, z) and Q(A,T, zq, z) are 
displayed for the initial condition A = 2 (Deuterium), zq = and different values of T, 
corresponding to energy E indicated in the figure. A(z) is evolved to Aq = 56 (Iron) and 
z g is determined. Then T g = Q(A,T , zq, z g ) is calculated as shown by the end points in 
the right panel of Fig. [U The trajectories for other values of A, e.g. A = 9, 14, 40, 
are very similar. The characteristic feature for all of them is the explosive regime in A{z) 
evolution at the end of the trajectory. We will discuss it in subsection 12.41 

In the calculations of the evolution trajectories we always include two additional con- 
ditions: maximum energy of acceleration E^ x (in terms of T^ x ) and nucleus stabil- 
ity restriction. The first restriction is imposed as follows. With z g determined from 
A(A, r, Zg) = Aq we calculate T g = G(A, T, z g ). If T g > this trajectory is forbidden, 
and the contribution to flux along this trajectory is set to zero. The second compulsory 
restriction, caused by nucleus stability, is considered in the next subsection. 

From the point of view of the final result of our work as calculation of UHE nuclei spec- 
tra, the evolution trajectories are needed to different extent. For trajectory method they 
are the most essential component of calculations. For the combined method the trajecto- 
ries are needed for calculation of two quantities, z g and T g . The former is calculated from 
A(T, A, z g ) = Aq (left panel of Fig. d]). The latter is calculated as T g = Q(T,A,z g ) (the 
end-points of trajectories in the right panel of Fig. 0]). In the coupled kinetic equations 
the A(T, A, z g ) trajectory formally is not included, but it presents implicitly, providing 
the evolution from Aq to A. 



2.3 Nucleus stability restriction 

A nucleus A with Lorentz factor V < at z = is stable, i.e. it is not photo-disintegrated 
during the Hubble time. This is easy to understand from the plots in Fig. [5j 




Figure 5: Stability condition for A = 40 and A = 56 nuclei. The nucleus A is stable if the 
decay rate r^ 1 is less than the rate of diminishing the Lorentz factor r" 1 + H , where r ee is the 
characteristic time of the pair-production energy losses (see Eq. [8]). The critical Lorentz factor 
given by Eq. (jSJ) coincides approximately with the intersection of r^ 1 and r" 1 curves, as shown 
in the figure. 

The critical Lorentz factor is defined by Eq. ([8]) equating the photo-disintegration 
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lifetime ta = (dA/dt)~ l and the characteristic time rr for changing of the Lorentz factor. 
Consider as an example the forward-in-time evolution of 40 Ca nucleus with the Lorentz 
factor r < T c . From Fig. [5] it is clear that such nucleus is stable: decreasing of the 
Lorentz factor goes faster than changing of A and the nucleus slows down, with a further 
increase of ta(T). Thus, nuclei A — 1 with T < cannot be produced by decay the 
A — > (A — 1) + N at z = 0. However, they can be produced at larger z. 

Consider now the backward-time evolution of the nucleus (A — 1) with Lorentz fac- 
tor r < . The Lorentz factor of this nucleus increases with z quite fast as T(z) = 
Q(A — 1, T, Zq, z), while A — 1 remains practically unchanged, because r^ 1 <C T~ ail . Simul- 
taneously rjr(z) decreases with increasing z, because the intersection point of t^ 1 and 
r pair * s shifted t° lower energies by a factor (1 + z). As a result at some critical redshift 
z c (T), determined by the condition ta(T,z c ) = rr(r,z c ) the considered nucleus becomes 
unstable. 

In case of nucleus A produced in (A + 1) — > A + N decay the equation for critical 
redshift z c becomes 

g(A,r,z ,z c ) = rf +1 (z c ), (11) 

where zq = 0. Therefore, to produce the secondary nucleus A with Lorentz factor V < 
r^ +1 at zq = 0, the parent nucleus ^4 + 1 must decay at z > z c determined by the Eq. 
(jlip . Thus, z c (T) can be considered as the minimum red-shift z m i n for the production of 
nucleus A with Lorentz factor T < 

A trajectory for A-nucleus can be blocked not only by A + 1 nucleus, but also by 
A' > A + 1, if A'-nucleus has an anomalously high T c . But in this case the nucleus A' — 1 
increases its Lorentz factor with fixed A' — 1 until it reaches r c . Numerical calculations 
show that stability condition for the considered case puts no restrictions to the calculated 
trajectories. 

2.4 Explosive regime in A(z) evolution 

The most striking and important feature of A{z) evolution, seen in the left panel of Fig. [H 
is its explosive behaviour at the end of the evolution. It is explained by the short lifetime 
ta of nucleus A relative to the transition A — > (A + 1). This lifetime is connected with 
f3° A , calculated for z = and plotted in Fig. [2J see also Eq. ([3]), as 

ta(T, z) = (1 + zr'A- 1 /^ + z)T\. (12) 

This equation demonstrates that ta(T,z) is short and diminishes with z. Taking into 
account that cr^is in Eq. ([7j) is approximately proportional to A one can understand that 
A(t) increases exponentially with time. The Lorentz factor at this stage remains constant 
because the ratio ta/ty is less than 100 (see Figs. [T]and[2|). 

At large Lorentz factors r > T c the explosive stage starts at very small z (see the 
curves 10 19 and 10 20 in the left panel of Fig. and the Lorentz factor during this stage 
remains constant (see right panel of the figure). 

Let us consider now a nucleus A at zq = with Ta *C This case is shown for 

A=2 in Fig. H (left panel) for T = 1.5 x 10 8 , i.e. E = 3x 10 17 eV. From z = to z ~ 0.5 
a nucleus evolves not changing A and changing T only due to adiabatic energy losses. The 
critical moment z c occurs when Lorentz factor reaches T c for 3 He according to Eq. (jlip . 
Then photo-disintegration starts and almost immediately enters the exponential regime 
with ta given by Eq. (|12p . which in terms of redshifts is an explosive regime. 
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Figure 6: The explosive regimes in the zoomed scale. Left panel: The explosive regime for 
nuclei with different A and with initial energy 3 x 10 17 eV at zq = 0, shown in a natural time 
scale. The numbers give A. The duration of explosive regime is seen to be short. Right panel: 
The explosive part of v4(,z)-evolution for A = 24 and E = 3 x 10 17 eV at z — 0. It starts from 
A—t(A+l) transition at za = 5.58 and finishes at z g = 5.59 reaching Aq = 56. 



The typical cases of explosive trajectories with Fa <C r^ +1 at zq = are shown in Fig. 
[U In the left panel the explosive regimes for A = 2, 14 and 40 are shown in a natural 
time scale to demonstrate the short time of evolution. In the right panel the evolution of 
nucleus A = 24 is shown for E = 3 x 10 17 at z = 0. Until z = 5.58 A is not changed, and 
then during Az = 0.01 A jumps to Aq = 56. It illustrates a typical case when the redshift 
of A (A + 1) transition approximately coincides with the redshift of Aq appearance, 
za ~ z g . In this example za = 5.58 and z g = 5.59. Considering then the part of trajectory 
A — > Aq, one obtains T g = Q(A, Ta, za, z g ) w Ta, due to za ~ z g . The relations 

z A ~ z g , T A w r ff . (13) 

are important feature of the explosive trajectories. 

For the large Lorentz factors T > T c , Eq. ([13]) is naturally fulfilled, too. 



2.5 Fluxes of nuclei and protons in the trajectory method 

In principle trajectory method is similar to Monte Carlo (MC) simulation: both are fol- 
lowing trajectory of a nucleus propagating through CMB, calculating the flux of secondary 
nuclei with Lorentz factor T from the flux of primary nuclei Aq with Lorentz factor T g . 
As advantage MC includes fluctuations in interaction of nucleus with CMB and does not 
need to introduce the limits of integration. 

In the trajectory method the space density n a (T,ZQ), of the components a = Aq,A,p 
is found from the conservation of the number of particles as n a (T)dT = J dt'Q a (V ,t')dT' , 
where Q a (T',t') is the generation rate of these particles (see appendix [A]). In a way 
convenient for calculations this formula has the form 



n a (r,2:o) = / dz' 



dt 

dz 



Q a (T',z')^, (14) 



where zq = 0, Q a (T',z') is the generation rate of a particles given in appendix lAl and 
dT'/dT is given in appendix [Bj The limits of integration follows. 



For primary nuclei (a = Aq) z m \ a = and for z max we have three options: z{ja ), 
i.e. redshift corresponding to lifetime ta (T), z c (T) as given by Eq. (fTTjh and z 9 (Aq + 1). 
These three cases give somewhat different results. 

For secondary nuclei (a = A) z m ; n = for T > T c , and is z c (T) for T < T c . The 
maximum limit is given by za+i(T). 

For the secondary protons produced by (A + 1) nuclei the limits are the same as for 
secondary A-nuclei, since the production of both is identical (see appendixes lAl and |D|) . 

As numerical calculations show, the trajectory method gives the worst accuracy for 
flux calculations, especially for secondary nuclei, and this is most probably connected with 
rather rough estimate of the integration limits, which affect strongly the results. 



3 Kinetic equation combined with trajectory cal- 
culations 



In this section we calculate the fluxes (space densities) of primary nuclei, secondary nuclei 
and protons using the kinetic equation and trajectory calculations. The latter are used 
in a limited way, only for calculation of generation rate of nuclei and protons, using the 
number of particles conservation (see appendix |A]) . In the case of secondary nuclei, the 
trajectories are used only to determine T g and z g associated to the primary nucleus Aq, 
these calculations may be greatly simplified taking into account that in all practical cases 
the trajectories are explosive. 

We start with the calculations of the secondary-nuclei flux as the most general and 
technically most interesting part of our calculations. 

3.1 Secondary nuclei 

We consider an expanding universe homogeneously filled by the sources of accelerated 
primary nuclei Aq with a generation rate per unit of comoving volume Qa (F, z) given by 



where j g > 2 is the generation index, is the nucleon mass and jCq is the source 
emissivity. i.e. the energy generated per unit of comoving volume and per unit time at 
z = 0. In Eq. (fT5j) r m i n ~ 1 is assumed. In all calculations below we assume also a 
maximum energy of acceleration E^ x = Zq x 10 21 eV (or = (Zq/Aq) x 10 12 ) with 

the condition Q^ (r g ) = at T g > T^ x . 

The trajectory A(z) is calculated as a continuous quantity only to determine the 
redshifts of (A + 1) — > A and A — > (A — 1) transitions. Between these values we assume 
A = const calculating the evolution T(z) with fixed A. Thus, we assume that both of the 
above-mentioned transitions are instantaneous photo-disintegration, e.g. 7cmb + (^.+1) — > 
Since a recoil momentum in these processes is negligibly small, one has approximate 
equality of Lorentz factors of all three particles 



The generation rate of secondary nuclei (and secondary nucleons) is derived in ap- 
pendix [A] using the trajectory calculation combined with the conservation of the number 
of particles 




(15) 



(16) 
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Qa(Ta,z) = Q Ao (T g ,z g )^- ^ (17) 

i -f- Zg ai a 

with dTg/dTA given in appendix iBl 

The kinetic equation for the comoving space density of secondary nuclei A, ua(Ta, t), 
reads 

dn A (T A ,t) d n A (T A ,t) n , , , , 

w, ~r- M r A, t)nA{T A , t)\ + — — = Q A (TA,t), (18) 

dt dT A T A (T A ,t) 

where 6 a = — dT /dt = (/3 pa i r + 13^)^ a is the rate of the Lorentz-factor loss and Qa(^a, t) 
is the generation rate given by Eq. (|17p . the nucleus A lifetime is 



1~A ( r A,z) 



-dA/dt 



2vr 2 r 2 



^o(A) 



dea dis {e, A)v{e)e 



In I 1 — exp[- 



2TT 



(19) 



where v{e) is the average multiplicity of the emitted nucleons. 

The solution of the kinetic equation (|18p is found in appendix ICl and can be presented 



as 



n A (T) 



dzA 



dt A 



dz A 



QA(TA,ZA)^-e~ 



(20) 



where ua(T) is given for zq = 0, za is the redshift of A-nuclei production in ^4+1 — > A+N 
decay, Ta = Q(A, T, zo, za) is the Lorentz factor of the nucleus A at the moment of 
production, calculated using A = const. The quantity t](Ta,z) is given by 



dz 



dt' 
~dz' 



1 



t a (T a ,z') 



(21) 



The physical meaning of the factor exp(— rf) becomes clear from Eq. (j76|) of appendix 
[Cl Assuming ta = const and introducing as variable a propagation time t = to — t' , one 
obtains this factor as exp(— t/r A )i which gives the survival probability for the nucleus A 
during the propagation time t. 

The factor exp(— rf) provides an upper limit of integration in Eq. ()20p . Without it the 
trajectory of integration in kinetic equation Ta{z) extends to any large z with the same 
A = const. Considering this trajectory as the one on which A changes with z too, we can 
put the marks za+i and zq, which corresponds A(z) reaching (A + l) and Aq, respectively. 
From Eq. (]2ip one can see that at these redshifts the suppression factor exp(— rf) is equal 
to exp(— 1) and exp[— {Aq — A)]. We shall use za+i and za as trial upper limits in Eq. 
(I20p for comparison with the basic case when the upper limit is set to infinity. 

The calculation of the lower limit of integration in Eq. (|20p is more complicated. The 
lower limit z m i n {T) is provided by vanishing of the generation rate Qa{z) at z < z m in- As 
has been discussed in section [2731 at z=0 there is a critical Lorentz factor below which 
the photo-disintegration of the parent nucleus A+l is absent (see Fig. [5]). With increasing 
z this critical Lorentz factor F^ +1 {z) slightly diminishes because both the pair production 
and photo-disintegration curves are shifted in Fig. [5] by factor (1 + z) to lower energies. 
If Lorentz factor of a considered A- nucleus Ta < the transition A — >• (.A + l) at 

z = is forbidden and only T increases with z {t^ 1 in Fig. [5] is very small). At z c when 
the Lorentz factor of A-nucleus G{A, T,0, z c ) reaches T^ +1 {z c ) photo-disintegration starts 
(see Eq. dTJ. Thus, for T < T^ +1 z min {T) = z c {T), and for T > T c z min {T) = 0. 
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Figure 7: Redshift limits in Eq. (|20|) for nuclei with small mass number A. The values 
and Za correspond to redshifts when the running A(z) reaches A + 1 and A , respectively. 




This effect of z m ; n -appearance can be explained in the evolution-trajectory approach 
as follows. The production rate of ^-nuclei, Qa(z), is caused by the decay of A + 1 
nuclei. However, at all z for which Ta+i(z) < T^ +l (z) a nucleus ^4 + 1 does not decay 
and A + 1 = const, as we observe indeed in the calculation of the trajectories forward in 
time, and hence Qa{ z ) = at z < z c . 

The lower limit z m i n (r) and two trial upper limits are presented in Figs. [7]and[8l The 
upper limit za+i gives the redshift when A formally evolves to A + The upper limit za 
corresponds to the evolution of A to Aq. The actual upper limit in this case is provided 
by the factor e _T? . One may recognize in these figures the regions in the energy-redshift 
plane which contribute significantly to the flux of the secondary nuclei. In particular, at 
low energies the interval z max — z m \ n tends to zero and the flux vanishes. At the highest 
energies E ^> 3 x 10^ cV, the interval is widening, but the absolute values of ^min 

and 

z max become vanishingly small, in the range of 10 -3 , and thus the fluxes are strongly 
suppressed. 

We will now come over to calculation of fluxes in the form of particle densities re^r, z). 
For this we calculate the generation rates Qa(^a, za) using Eq. (fT7|) and put it into Eq. 
(|20p . The left panel in Fig. [9] explains the procedure of integration. The fluxes are shown 
in Fig. [10] in terms of energy (left panel) and Lorentz factor (right panel). 




Figure 9: A sketch of T(z) evolution in kinetic equation approach. Left panel: The thick red line 
gives the integration trajectory in Eq. f|20|) described by T^(r, z , z) with fixed A = const. The 
auxiliary A-evolution trajectories (A + 1) — > A are needed for calculation of generation rate 
of A-nuclei Qa(^a, za)- In practice all these trajectories are explosive, providing za ~ z g and 
T g . A point (za^a) represents the running variables in the integral given by Eq. (120]) . 
The lower limit of integration is given by z^ m , while za+i and za are the trial upper limits 
(see text). Right panel related to section [372] The thick red line gives the integration trajectory 
in Eq. (I2rj|) for secondary protons. The trajectory describes the Lorentz-factor evolution of a 
proton r^(r p , z). Proton is produced at a running point (z, T) of a proton trajectory in a decay 
(A+1) — > A + N, and we use in the figure the variables of the brother nucleus A: za = z 
and Ta = T. The A-evolution sub-trajectories A — > Aq are needed for calculation of generation 
rate (Ta, za) given by Eq. f l27|) . All these trajectories are explosive, providing za ~ z g and 
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One should realize the following feature of calculations: the time scale of A evolution 
from (A + 1) to Aq needed for calculation of the rate Qa(^a, z a), and the time scale of 
the kinetic-equation trajectory Ta(z), which contribute the integral in Eq. (|20p. is the 
same and very short ~ T4. In case of ^4-evolution this is the explosive part of trajectory, 
in the case of kinetic-equation trajectory it is regulated by exp(— t/r^). For low energies 
E < 10 18 eV the lower limit of integration z m \ n = z c ~ 2 — 3, and A evolution at this 
z occurs in explosive regime with short evolution time given by Eq. (|12p . The region of 
integration in Eq. (|20|) is controlled by the same ta(I\ z). As it is easy to understand 
the situation is the same at V > T c . In both cases we can use an explosive trajectory for 




Figure 10: Fluxes of secondary nuclei, as function of energy (left panel) and of the Lorentz factor 
(right panel). The numbers on the curves show A. The small filled circles present calculations 
in the explosive trajectory approximation (see the text). 



calculation of generation rate, which implies z g ~ za and T g ~ Ta, and thus Eq. (|20p 
results in 

Mr) = ^r ,. t " Z * i f 7"^^\ (22) 

m N A J Zrnin (r) (1 + za)H(za) dT 

where Ta = Ta(T,zq,za) with A = const and zq = 0. 

In Fig. [10] the fluxes calculated in this explosive trajectory approximation are shown 
by small filled circles. One may observe, as expected, the good agreement with the exact 
calculations. 

The basic parameter which determines the spectra of nuclei A is given by the critical 
Lorentz factor T^ +1 at which the photo-disintegration process of the parent nucleus (A+l, 
or Beryllium in the case of Helium) is allowed at z = 0. One may expect that the peaks 
seen in the spectra E a tia(E) set at the critical Lorentz factors and the right panel of Fig. 
[10] confirms this expectation with a good precision (see Table [TJ of section [2TTT) . 

The calculated spectra at low T exhibit the A hierarchy: the higher A the larger the 
flux. The exceptional case of low Be flux is explained by large ta, which determines the 
upper limit in Eq. (j20|) . In effect, the upper limit za+i is below z m \ n (see upper-right 
panel in Fig. [7]), and factor e~ v in Eq. (|20]l operates exclusively at large n suppressing 
strongly the Beryllium flux. 

The explosive-trajectory approximation allows us to obtain the exact formula for the 
asymptotic high-energy flux and to compare it with calculations presented in Fig. [TOl 
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In the asymptotic regime one can use in Eq. (|22p Ta ~ r (Lorentz factor is not 
changing in the pure photo-disintegration process), z m \ n = and z g -C 1. Then after 
simple calculations we obtain 

n A {T) = ( lg - 2)^-T~^r A {T). (23) 

The asymptotic formula (|23h follows also directly from the kinetic equation (I18p if one 
neglects there dn A /dt and <9(&a"-a) / 'dT a- Eq. ([23]) predicts the asymptotic ratio as 

n Al ■ n A2 ■ riA 3 - = ta x ■ ta 2 ■ ta 3 -- (24) 

Numerical calculations confirm that asymptotic relations (|23ft and (|24ft are fulfilled for 
the curves in Fig. [TUJ 

For all secondary nuclei at the highest energies the contributing universe is vanishingly 
small, when the relevant redshifts are of the order of 10 -3 or less, as it follows from curves 
za+i{E) in Figs. [7] and [8] The fluxes in Fig. [10] are non- vanishing only because in the 
present paper we consider a homogeneous distribution of sources for primary nuclei Aq 
with constant density at any redshift, while in a realistic situation the presence of UHECR 
sources at these low redshifts (z < 1CP 3 ) is unlikely. 

The energy spectrum at high energies follows the corresponding energy losses, with a 
recovery at the highest energy due to an increase of the nucleus lifetime ta (i-e. a decrease 
in the photo-disintegration energy losses, see Figs. Q]and[2]). 



3.2 Secondary protons 

In this section we discuss the flux of secondary nucleons produced in the photo-disintegration 
of the heavier nuclei. We do not need to distinguish neutrons and protons, because at 
the characteristic lengths involved here neutrons decay fast at all energies of interest. For 
this reason we will often refer to protons, instead of nucleons. 

Production of secondary protons accompanies the production of secondary nuclei, but 
in contrast to them, protons are not destroyed and can arrive from any redshift, being 
suppressed only by high energy of generation. 

Following the same approach as for secondary nuclei one can write a kinetic equation 
that describes the propagation of protons: 

_ JL [ 6p (r p , t )n p (T p , t)] = Q P (F P , t) (25) 

where n p is the secondary-proton density, b p = —dT/dt describes the loss of proton Lorentz 
factor (energy loss) due to adiabatic energy losses, pair-production and photo-pion pro- 
duction on the CMB radiation, and Q p is the generation rate for secondary protons, 
produced by photo-disintegration of secondary or primary nuclei. The secondary nucle- 
ons N which accompany production of ^4-nuclei in the process A + 1 — > A + N have the 
same Lorentz factor and generation rate as A-nuclei. We refer to these nucleons as U A- 
associating protons" with notation for their space density, and use for their generation 
rate Qp(T,t) = QA(T,t). The solution of Eq. (j25|) for ^-associating protons similarly to 
[52] . reads: 



n 





dt' 


/.. d ~ 


' ^min 


dz' 



1 



where T'(z') = G p (T p , 0, z') and dT'/dT p is taken along the proton trajectory ( T p g (T p ,z) 
line in the right panel of Fig. [U] ) . 

As in Eq. (|17p we have from the conservation of particles number: 

where the redshift of A + 1 decay za = 2' and Ta = T' , index A at derivative dT g /dT' indi- 
cates the A- variable subtrajectory A — > Aq for which derivative is given by Eq. (|68p . and 
Qa is th e generation (acceleration) rate of primaries Aq at redshift z g with Lorentz factor 
r 9 , given by Eq. (fT5|) . The generation redshift z g is determined from A(A, T', z' , z g ) = Aq 
and the generation Lorentz factor is calculated as T g = Q(A, T', z', z g ). Finally, we obtain 
for rip (T p , zq) at zq = 0: 

v a (t v _ i 9 - 2 , dt> i + z> /dv g \ ( dr \ 

the integration in Eq. (|28|) goes along the proton trajectory T p (T p ,z), i.e G p (T p , 0, z), as 
sketched in the right panel of Fig. [9l A proton can be produced at any running point 
of this trajectory (za,Ta), with za = z' and Ta = T', in a decay (A + 1) — >■ ^4 + AT. 
The generation rate at such point Q p (za,Ta) is calculated using Eq. (|27|) along the A- 
evolution subtrajectory A Aq. For the different subtrajectories we have different values 
of T g for the primary nucleus Aq and different redshifts of generation z g . 

The minimum redshift in Eq. (|28[) is calculated as explained in sections 12.31 and 13.11 If 
the proton Lorentz factor T p < T^ +1 (z = 0), the minimum redshift corresponds to z c at 
which the proton Lorentz factor T^(z c ) reaches the critical Lorentz factor of the (A + 1)- 
nucleus T^ +1 (z c ). Note that z c = z m j n for protons are higher than for secondary nuclei 
(compare Fig. [TT]with Figs. [7]and[8|). It happens because /3^ ir (T) is Z 2 /A w A/4 of that 
for proton (see Eq. [6]). In the case T p > T c , z m ; n = 0. The maximum upper limit z max 
corresponds to A(z) subtrajectory on which the Lorentz factor of the ^o-nucleus reaches 
T™^ x (see Fig. [9]). These limits are plotted in Fig. [TTJ with z max shown by solid red curve. 
The other upper limit shown by small filled circles is valid for the instantaneous method 
(see below). Here it is important to note that both limits coincide with great accuracy. 

The total space density of secondary protons is given by the summation over the 
primary Aq and all the secondaries A in Eq. ()28f) : 

n p (T p ) = <(r P ) • (29) 
A<A 

Before discussing the numerical results for the calculated flux, we present an alternative 
and more simple method of calculation, which, being less precise, can be considered as a 
test at high energies for the calculations described above. 

This approach will be referred to as "primary-nucleus instantaneous decay". It is based 
on the assumption that at the moment of production z g the primary nucleus is instanta- 
neously photo-disintegrated to nucleons. At large Lorentz factors this assumption is 
well justified because the nucleus lifetime TA (T g ) is much shorter than all other relevant 
time scales of the problem. We will demonstrate below why this approach works well 
(with some exceptions) at lower Lorentz factors r < T c and T <ti T c . To distinguish the 
previous more precise approach from the instantaneous one, we will refer to the former as 
intermediate A-decay method. 
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Figure 11: Upper and lower redshift limits in Eq. (|28l) . The lower limits are shown for three 
proton-associating nuclei with A=4, 9 and 56 and labelled by these numbers. The highest 
limit for Helium and the lowest for Beryllium are explained by the highest and lowest critical 
Lorentz factors T c for these nuclei. The lower limits for A= 14, 24, 35 and 40 (not shown here) 
are grouped between A=4 and 56. The upper limit (see text) is given by full red curve. The 
small filled circles show the upper limit for the instantaneous approach. 




The proton space density in the instantaneous approach can be written as 



n p ast (T p ) = A Q / dz g 



dtg 



dZg 



dT p 

Q Ao (TP,z g )-^-, (30) 



v 



where all notations are as before and Qa are given by Eq. ([15]) . 

Since Qa oc 1/Aq, the factor Aq disappears from Eq. (f30j) and n p nst (T p ) can depend 
on Aq only through the limits of integration. 

After a simple rearrangement of Eq. (|30p one obtains 

,_^0 7 9 - 2 ff" ^ [TP(F p , z g )]-T dFP 



< st (r P ) = / dz 9 9 *' y " (3i) 



Eq. ()3ip is very simple for calculations, because it involves only the proton trajectory in 
r — z plane (subtrajectories A — > Aq in the right panel of Fig. [9] shrink to a point). The 
limits of integration are the same as those shown in Fig. [TTJ The basic observation involved 
in the calculation of these limits is that the Aq nucleus is now located on the proton 
trajectory TP(T p ,z) in Fig. Therefore, z™ m is defined by the equation TP(T p , z m - m ) = 
r^°(-Zmm)) which is the same as for nuclei A, considered now as Aq in Fig. [TTJ The upper 
limit is defined by the condition r^(r p , z max ) = F™^ x and it is shown in Fig. [IT] by small 
filled circles, which coincides with upper limit in the intermediate A method. We see, 
thus, that the upper limit z max does not depend on ^o- The lower limit at T p > is 
Zmin = 0. At T p < T c the lower limit depends on Aq, but weakly, with the exceptional 
case of Beryllium. 

Thus, at high energy when z m i n = the flux has a universal form independent of the 
primary nucleus, i.e. of Aq. It is natural to expect that the high energy regime starts 
from T p > when z m i n = 0, and the calculations below confirm this expectation. 

The high-energy regime is of particular interest as a test of the intermediate A- decay 
method. In this case the proton flux in instantaneous approximation n p nst must coincide 
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Figure 12: Spectra of secondary protons computed by the A- intermediate method (continuous 
red curves) and in the instantaneous approximation (black dotted curves) for different A Q . The 
generation spectral index is 7 S = 2.3 in all cases. 




Figure 13: The same as in Fig. [T2]for heavy nuclei. 
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with the intermediate A method, when the flux is given by J2 n p > i- e - by Eq. (|29p . 
It follows also from Figs. [T] and [21 which show that at T > T c the photo-disintegration 
rate strongly dominates over the rate of the Lorentz-factor loss. Figs. [12] and [13] reli- 
ably confirm the agreement of both calculations at T > r c , demonstrating thus that the 
intermediate A method has passed this test. 

We address now the question of the agreement between A-intermediate and instan- 
taneous methods at T < T c and r < r c . The explanation follows from the explosive 
character of A(z) trajectories for A — > Aq evolution. Using the explosive-trajectory result 
Z A ~ Zg an d ss in Eq. (|28|) . one obtains after simple calculations 

<{Tp) = f " dZa - l ^£L ■ g (32) 



The total flux of protons is obtained by summation in Eq. (|32p over all A. Since limits of 
integration depend weakly on A, it means multiplication of Eq. (j32f) to Aq, which results 
in instantaneous flux given by Eq. (|3ip . 

We shall discuss now the calculated spectra presented in Figs. [T2l and [T3l 

One may observe that agreement between the exact (^4-intermediate) and instanta- 
neous method of calculations is precise at high energies (r > r c ) and approximate at low 
energies ( (T < T c ) and (T <C r c )), as it is expected. 

For small A the excellent agreement at low energies is observed for 4 He. This is an 
exceptional case: the energy threshold of photo-disintegration for 4 He is much higher 
than for its daughter nuclei: 3 He and D. Hence, when 4 He starts to photo-disintegrate, its 
daughter-nuclei disintegrate too. Based on this argument we expect that in the Helium 
case the spectra of secondary protons calculated by both methods should be identical 
at all Lorentz factors, including the low ones, and this expectation is confirmed by the 
upper-left panel of Fig. [12] 

We now come over to the "Beryllium excess" in the spectrum of the instantaneous 
decay, shown in the upper-right panel of Fig. [T2l The nature of this effect is simple: 
r^ e <C r^ e and since 4 He is born with Lorentz factor of 9 Be, the former can be stable and 
secondary protons are not produced. For V > T® e it gives a suppression of the secondary 
proton flux. Under the assumption of an instantaneous decay this effect is absent. 

Beryllium effect works for all Aq > 9 and since it affects only 8 protons from two 4 He 
nuclei, its influence on heavy nuclei is weaker, as one observes from Fig. [TBI 

The shape of the proton spectra in Figs. [12] and [13] can be naturally explained. They 
have two spectral breaks: the high energy steepening, which is the usual GZK cutoff at 
energy E ~ 5 x 10 19 eV, and the low-energy steepening which coincides with T c : ~ 4 x 10 9 
for 4 He, ~ 4 x 10^ for 9 Be etc). Below r c the adiabatic energy losses dominate and it 
explains the flat spectrum n p (E p ) there, i.e. steep E 3 n p (E) spectrum in Figs. [T21 and 1131 
It is different from oc because z m i n in Eq. ([28]) and Eq. ([3T]) depends on T, i.e. on 

E p . 

We complete this section with a remark about the universality of the secondary-proton 
spectrum. Using Eq. (|3ip for the instantaneous spectrum we proved its universality, i.e. 
independence of Aq. The intermediate A spectrum has the same property, as it was proved 
above. In other words the normalised spectrum calculated for one Aq is valid for any other 
Aq with the same normalization. 
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Figure 14: Spectrum of primary nuclei with various Aq as function of energy (left panel) and 
Lorentz factor (right panel). 



3.3 Primary nuclei 

We calculate now the diffuse flux of primary nuclei Aq, which arrive undestroyed from the 
sources. The comoving space density of these nuclei n Ao (T,t) is described by the kinetic 
equation 

dn Ao (T,t) d n Ao (T,t) 

di df ^ ^ ' ^ " t a (r t) = ^ ' ' ^ ' 

where Q Ao (T,t) is the generation rate of primary nuclei per unit time and comoving 
volume, and T Ao (T,t) is the photo-disintegration lifetime of the nucleus Aq calculated as 
[A p dis (A ,T, z)]- 1 with /3 dis given by Eq. 0. 

As found in appendix O the solution of the kinetic equation (|33|) reads 



n Ao (T) = / dz 



oo 







dtg 



dZg 



QA {T a ,z g )^e-^ T ^ (34) 



where n Ao (T) is the space density of nuclei Aq at z = 0, T g = G(Aq,T,zq = 0,z g ) is the 
Lorentz factor of nuclei Aq at the moment of generation z g calculated by the evolution 
trajectory with fixed Aq, which starts from r at zq = 0. The ratio dT g /dT is also 
calculated with fixed Aq. The quantity r\ takes into account the photo-disintegration of 
the propagating nucleus: 

. . f to dt . , 

V{T g ,z g )= . (35) 

Jt(Zg) TA (1 (t),t) 



The lower limit of integration in Eq. (|34p reflects the assumption of an homogeneous 
distribution of the sources. The upper limit is imposed by the factor exp(— rj) accompanied 
by the condition of a maximum acceleration Lorentz factor: Q Ao (T g ) = 0, if T g > T^ x . 

In Fig. [T4]we plot the flux of primary nuclei for various Aq, using the injection spectrum 
with 7 9 = 2.3. The steepening of the spectra in Fig. [T4] starts from the Lorentz factor 
where pair-production energy losses are equal to adiabatic energy losses. Comparison of 
Figs. [5] and [T4] confirms it. It is worth noting that this criterion for the beginning of the 
steepening was obtained first in [31] . 



4 Coupled kinetic equations (CKE) 



In this section we present a straightforward method of calculation for UHE nuclei fluxes, 
based on the analytic solution of a complete set of coupled kinetic equations. The great 
advantage of this method is a simple expression for generation rates of secondary nuclei 
and protons, which do not include integration with independently determined lower and 
upper limits. We consider as before the expanding universe filled homogeneously by 
accelerated primary and secondary nuclei, including the secondary protons. As discussed 
in section 12.11 these secondary particles are produced mainly by the photo-disintegration 
process (A + 1) + 7 — >• A +p. In this approximation one can write the generation rate for 
any secondary nucleus A as a function only of the equilibrium distribution of the parent 
nucleus (A + l), namely: 

Q A (T, z) = — — , (36) 

TA+l(l ,Z) 

where t a (T, z) is given by Eq. (|12p and (|7|). Taking it into account one may write a set of 
coupled kinetic equations which describes the equilbrium of primary nuclei Aq with the 
products of their decay 

dn A{) (T,t) d n Ao (T,t) 

—at ^ [nAo(r '^ (r ' t)] + 7^) = QA ^ t] 

dn Ao -i(T,t) 8 (r +\ hA (r ^1 , "A,-i(r,*) n Ao {T,t) 

! (37) 

9n A (T,t) d n A (T,t) n A+1 (T,t) 

7- ^ [n A {T, t)b A {T, t) + — — — = — — 

ot dT T A (r,t) TA+i{T,t) 

The solution of the kinetic equation for primaries Aq is given by 

P Zmax r\ ffX^^ 

- (r ' 2 » = / (irki e * (r '^^' (38) 

with 

f z ' dz" 1 

^ T '^ = j z (iTW)^' (39) 

and the solution for an arbitrary secondary nuclei A is 

n A{T z) - f ZmaX d A n A+l {T^z') dV {v , y) 

A{ ' j " J z (1 + z')H(z') r A+1 (r, z>) dT ' m 

with 

f z ' dz" 1 
r] A (T',z f )= - — — . (41) 
J z (1 + z")H(z") t a (T , z") 

In each of these equations we assume A = const and the Lorentz factor changing due to 
adiabatic and pair-production energy losses. Accordingly, the ratio dT'/dT is given by 



dT' 1 + z' 

W = TT7 exp 



Z 2 r z ' {l + z") 2 dz" (db p (ty 



A J z H(z") \ dT 



r=(i+ 2 ")r" 



(42) 



no 



where %(Y) = dT/dt is the Lorentz-factor loss per unit time for protons at z = due to 
pair production. 

The procedure to solve the system of equations (|37p consists in finding the analytic 
solution riA (r,t) for the first equation in the form (|38|) . and putting it into the second 
equation. Taking the solution of this second equation, given by Eq. ([4"0]) with A = Aq — 1, 
one finds riA Q -i(Y,t). Continuing this procedure one obtains the space density (fluxes) 
for all secondary nuclei nA(Y,t). 

Note, that in this method there is no problem with the limits of integration: the lower 
limit is z by definition and the upper limit is regulated by the factor e _r? , while z max here 
plays only a symbolic role. The generation rate Qa(Y, z) given by Eq. ([36]) tends to zero 
at small T due to factor ta+\ — ^ oo at small Lorentz factors. 

For the calculations of the secondary protons spectrum we use the ^4-associating pro- 
ton approach developed in section 13.21 The secondary nucleons N which accompany 
production of ^4-nuclei in the process {A + l) — > A + N have the same generation rate 



The kinetic equation for the secondary protons is given by Eq. (|25p with a solution 
described by Eq. (|26p . Using the generation rate in the form of Eq. (|43p this solution 



where Y'(z') = G P (T, z, z') is the proton trajectory, and Eq. (|44p is valid for any arbitrary 
initial z including z = 0. In fact, the actual lower limit of integration z m j n is automatically 
provided by ta+i(Y',z'), and this gives a great advantage of the CKE method in com- 
parison with the combined method, where z m ; n is calculated independently, using some 
particular features of the trajectory evolution. As will be seen by example of the sec- 
ondary A-nuclei, the brothers of the secondary protons, the structure of the lower-limit 
cutoff exposed by Eq. (|52p is quite different from the sharp cutoff at z m [ n introduced in 
the combined method. 

The upper limit of integration z m3iX in Eq. (|44p is quite different from the secondary 
nuclei, where it is imposed by the life-time ta in the form of e~ riA . In the case of secondary 
protons the upper limit is provided by r max , which the brother ^4-nucleus (or parent A + l 
nucleus) is allowed to have. Since at these very large Lorentz factors the approximation 
of explosive trajectories is fully justified, we can use r max for primary nucleus Aq and 
z max is given (in case of initial z = 0) by the equation for the proton trajectory r max = 

G P (r,o 

j z max ), where Y is the proton Lorentz factor at z — 0. 
As we discuss below, the fluxes obtained by the CKE method for secondary nuclei 
and protons differ, most notably at low energies, from that in the combined method (see 
section [3]). This difference is due to the lower limit of integration in the Eq. (|40p and 
d33J . In the combined method it is sharply fixed at z m i n as discussed in section 13.11 while 
in the CKE method the integrand is automatically suppressed by the term ta+i(Y',z') 
(see discussion below). 

The fluxes of secondary nuclei and protons are displayed in Fig. [T5l 
We will discuss now the solutions obtained for secondary nuclei and protons. In par- 
ticular we will study these solutions in their analytic form to compare them with the 
combined method. 



Q£(Y,z) = Q a (Y,z) 



n A +i(Y,z) 
t a+ i{Y,z) ' 



(43) 



reads 




(44) 
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Figure 15: Flux of secondary nuclei (the primary Aq = 56 included) (left panel) and of the 
secondary protons (right panel) as function of energy. The numbers on the curves show A. 
Ha is related with the total flux of protons produced by the primaries Aq = 56, the numbers 
indicate fluxes for each A-associating proton. The discussion of these results and comparison 
with the combined method is given in the next section. 



First we obtain the analytic solutions of Eqs. (|37h in the high-energy asymptotic 
regime, when all ta are very short. We re-write the solution (|38p in terms of the cosmo- 
logical time using for e _,? the expression of Eq. (j76l) (see appendix ICj) : 



,.dr' r /•* dt" 



n Ao (T,t) = J^ dt'Q Ao (T',t')—ex.p - 



r Ao (r",t" 



(45) 



Introducing the propagation time t prop = t — t' as integration variable in Eq. (|45p and 
taking QA dT' /dT out of the integral as a slowly variable quantity, one obtains after simple 
calculations 

n Ao (T,t) = Q Ao (r,t)T Ao (T,t), (46) 

and 

Q Ao -i(r,t) = n Ao (r,t)/T Ao (T,t) = Q Ao (T,t). (47) 
Repeating these calculations for A < Aq we obtain for the high-energy asymptotic regime 

Q A (T,t) = Q A +i(T,t) = - = Q Ao (T,t) (48) 

and 

n A (T,t) = Q A (T,t)T A (T,t) = Q Ao (T,t)T A (T,t), (49) 

i.e. the same high-energy regime that we obtained in Eqs. (j23j) and ()24j) in the combined 
method. In particular, the equality ()48|) is the same obtained in the combined method for 
the explosive regime, which is provided by short t a like in the case above. 

We will now come over to the general case valid also for low T. To obtain the gen- 
eral solution n A (T,t) in the analytic form we put into Eq. (|40p the value n A +i(V , z') 
determined from the preceding equation 

nA+l(r ,Z) ~L (i + z")H{z") t a+2 (f», z») dr 6 ■ (50) 



Repeating this procedure with increasing A until we reach Aq, we obtain the general 
expression with many-fold integral 



n A (F,z) 



_ rZmax p — VA -i(Ti,Zi) 

n / 



(1 + Zi)H(zi)TAo+i-i(Ti,Zi 



ZmaX dz' QAo{F l Z '\. dr %-wrv^ 



[l + z')H(z') dT 

(51) 

where dT' /dT is taken along the ^-variable trajectory. Using Eq. (fT5|) for Qa (T', z') and 
assuming that the external integration results in some functions f(T,z) (we calculated 
these functions numerically) we obtain 



"^•°> = I d ^L (TT7w?) £ ^ [r ' (r '- 4 ' 2 ' )r7l ^ <1 ' 3 ' 

(52) 

with /(z) rapidly vanishing when z — > and z — > z max . 

Comparison of this solution with n^r, 0) given in Eq. (|22p for the combined method 
shows similarity and differences. The crucial one consists in z m i n . While in the combined 
method this is a well determined function of T, in the CKE method the lower limit is 
given by different values of z weighted by the function f(z). 

Numerical calculations show a very good agreement between the two methods at the 
highest energies and for all primary nuclei Aq at all energies. For secondary protons the 
disagreement is moderate, but for the secondary nuclei there is a strong disagreement 
at low energies. We see the main reason for this disagreement in the z m i n (T) calcula- 
tions. In the combined model z m i n is calculated from reasonable physics connected with 
trajectories. In the CKE method the lower limit appears automatically, without any as- 
sumption, directly from the kinetic equations. We think that as far as numerical results 
are concerned, this last method should be trusted more. 



5 Discussion and Conclusions 

In this paper we have studied comparatively three analytic methods to determine the 
diffuse spectra of ultra-high energy nuclei propagating through background radiations: 
(i) trajectory method, (ii) kinetic-equation method combined with trajectory calculations 
(the combined method), and (Hi) coupled kinetic equations (CKE). 
We summarize first our general approach. 

We have calculated the diffuse spectra of primary and secondary nuclei, and secondary 
protons in the case of homogeneously distributed sources. According to the propagation 
theorem [52] in this case the spectra do not depend on the specific way of propagation, 
and thus one may use rectilinear propagation, or not thinking about space propagation 
at all, imaging all secondary nuclei and protons filling the space homogeneously together 
with primary nuclei. 

Evolution in time is the basic feature of our problem. A primary nucleus Aq is born due 
to acceleration and, interacting with the background radiation, give rise to the secondary 
nuclei and protons. 

With the help of two coupled differential equations for the propagation of a primary 
nucleus in the backward time, i.e. increasing redshift, we found the evolution trajectories 
A{z) = A(A, r, zq, z) and T(z) = G(A, T, zq, z), where the first three arguments fix the 
initial conditions. In the evolution equations the atomic number A changes continuously 
with z, but we use the continuous trajectory A{z) only to determine the integer A positions 



on the trajectory: A, A + 1 etc. At these points we assume an instantaneous photo- 
disintegration, e.g. (A + 1) + 7cmb — > A + N. The Lorentz factor between the points 
of photo-disintegration is calculated assuming A = const. The trajectories are basically 
needed to compute the generation parameters z g and T g , at which A(z) reaches Aq. 
Apropos, the formal solution of the coupled trajectory equations ([9]) with continuous A 
results in the same T g and z g as in our basic method. 

The important feature of the calculated A{z) trajectories is their explosive behavior 
for all initial energies (as an example see the left panel of Fig. [4] ). At small T < T c a 
nucleus propagates with A = const until large red-shifts z ~ 1 and then 'explodes' to Aq 
on a scale of photo-disintegration lifetime r a ~ 10 5 — 10 6 yr, see Eq. (|12p . It results in 
Z A ~ z g anci ~ F 9 , where index A marks the values at the threshold of an explosion. 
Therefore, the generation rate for A nuclei approximately equals to that of primary nuclei 
A) : Qa(Fa, z a) ~ QA {^g-, z g) with z g za and T g « T^. However, the main attention 
in the paper is given to the exact calculation of Qa(^a, z a) from the number of particles 
conservation. 

We describe first the combined method. 

The space density of each species a = A, p, Aq and their energy spectra are calculated 
using the kinetic equation 

dnai l a,t) ~ 4r Mr a ,t)n a (T a ,t)} + ^4 = ^ 
ot dT a r a (T a ,t) 

where b a = —dT a /dt = (/3p air + P^^a is the rate of the Lorentz-factor loss, T a (T a ,z) is 
the photo-disintegration lifetime of particle a (for a = p, r a = oo) and Q a (T a ,z) is the 
generation rate of particles a. 

For secondary nuclei A and protons, produced in the photo-disintegration processes, 
the generation rate Q a {Pa^ z ) is found from the conservation of the number of particles 
along a trajectory, dN a = dN g , where dN a is the number of produced particles a and 
dN g is the number of generated primaries at acceleration. The relation between the 
generation rate of primaries Aq and of secondary A and p is given by Eq. (|6Q|) . where 
z g and Tg are calculated using the evolution trajectories. The solutions of the kinetic 
equations ([53]) are found analytically. 

An important physical quantity in our calculation is the critical Lorentz factor 
at epoch z = 0. It is determined by the equality of the rates for changing of T and 
A'- r di s 1 (Tc l ) = TpairC^c ) + -^Oj where Hq is the Hubble constant, describing here the 
adiabatic energy loss, and indices 'dis' and 'pair' are related to photo-disintegration and 
pair-production lifetimes, respectively. This relation provides a stability condition for the 
A nucleus at z = 0: V < F^. The critical Lorentz factor r c is a basic energy scale, which 
gives a key for understanding all processes considered here and which also explains some 
features in the calculated spectra. 

Therefore, the combined method includes three elements of calculations: the kinetic 
equation, which gives the density of particles n a (T), the generation rate Q a (T a ,z a ) and 
the limits of integration, especially z^ in . The first two elements are reliable components 
of calculations, with Q a reliably evaluated with the help of the number of particles con- 
servation; the third element z m \ n is less reliable since it involves additional consideration 
and assumptions. 

The trajectory method of n a (T) calculation is less reliable because of the great uncer- 
tainties in the limits of integration z m [ n and z max in Eq. (I14p . 

The coupled kinetic equations (CKE) method is based on the set of kinetic equations 
([37]) . where the first one describes the primary nucleus Aq, the second Aq — 1 etc down 



to the A of interest. The solution of each preceding equation nA'(T, z) gives the gen- 
eration rate for a successive nucleus A 1 as nA>+i(r, z)/ta>+i(T, z), and the long lifetime 
TA'+i(r, z) at small T and z automatically provides z m \ a in the solution of a kinetic equa- 
tion, while the lifetime ta* provides, through the term e~ v , the upper limit. The A(z) 
trajectories are formally not involved in this method, but actually they are. The term 
it<a /ta q describes the disappearance of nuclei Ao in the first equation of the Eqs. (i37|) . 
but the same term appears in the second equation as generation of (Aq — 1) nuclei etc. 
Equation (|5ip . which includes the product of terms with all A' from A to Aq, demon- 
strates it more clearly. Indeed, this equation implies the evolution of T(t) with A' = const 
regulated by exp(— t prop /TA')- Then A' disappears, giving rise to the production term in 
the successive kinetic equation for (A' — 1). However, there is an essential difference with 
the combined method. The product in Eq. (|5ip is time-ordered which, according to the 
theory of probabilities, means fluctuations taken into account. 

As we indicated above the weakness of the combined method consists in the procedure 
of defining z m i n , while in the coupled kinetic equations this limit is a natural feature of 
the kinetic equation itself. The lower limit in the combined method is introduced by the 
stability condition given by Eq. ljlip . which is quite natural. However, in fact one can 
use the other definitions of stability. As a plausible example we consider an alternative 
possibility. 

We used above in the combined method the cutoff in integration over z imposing the 
condition QA{T,t) = at t > t max in terms of the cosmological time t (see the solution 
in terms of t given by Eq. [431) . This condition implies that at t > t max (A + l)-nucleus 
is stable. Let us introduce now this condition through the probability of A + 1-decay: 



One can see that the factor introduced is the decay probability of the A + l nucleus. It can 
be easily seen in the case ta+i = const, when this factor is 1 — exp(— t prop /TA+i), where 
iprop = t — t g is the propagation time. The decay factor in Eq. ([55]) provides the soft 
cutoff of the generation function. However, this cutoff also differs from the one obtained 
in Eq. (pD . 

We compare now the numerical results in some details. 

The two methods give numerically the same spectra for all nuclei at high energies and 
for primary nuclei at all energies. In particular both methods give high-energy fluxes of 
the secondary nuclei described by identical equations (|49p and (123j) with a recovery of the 
spectra at the highest energies due to an increasing ta(X). 

The peaks in the spectra of the secondary nuclei located at T c in the combined method 
(see Fig. fTO] right panel) coincide with the positions of peaks in the CKE method (see 
Fig. [15] , left panel). This is a surprising result, because the critical Lorentz factor is not 
introduced in the CKE method. The fluxes above r c are the same, but they are much 
different below, which is mainly a result of the differences in z m i n . The peculiar behaviour 
of 4 He at low energies in Fig. \W\ (left panel), being confirmed by explosive trajectory 
calculation, is not seen in Fig. PT5l 

The fluxes of the secondary nuclei obey as a rule A-hierarchy: the heavier A, the larger 
the flux (see Figs. [TUl and [15]) . 

The main features of the secondary-proton spectra are the same in both methods. The 
^4-associating proton spectra are the same for all A, and the total spectra are universal, i.e. 
do not depend on Aq. In the combined model this is a direct consequence of the explosive 




(54) 



character of the trajectories. At low energies (E < 1 X 10 19 eV) in CKE the A-associating 
proton spectra are different for various A (see right panel in Fig. I15|) . This is related to 
a breaking of accuracy of the explosive trajectory approximation at low energies. 

The shape of the total proton spectrum is the same in both calculations. It has 
two spectral breaks. The high energy steepening is the usual GZK cutoff at energy 
B~5x 10 19 eV. The low -energy steepening coincides with T c . Below T c the adiabatic 
energy losses dominate and it explains the flat spectrum n p (E p ) there, i.e. the steep one 
in terms of E s n p (E) spectrum as plotted in all figures here. 

In the second paper of this series, the EBL background will be also included. In this 
case the high energy part of the nuclei spectra is determined by CMB, as considered in 
the present paper, while the low-energy part is affected mostly by Infrared, Visible and 
Ultra Violet radiations, which compose the EBL. 
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A Generation rates of primary nuclei, secondary 
nuclei and secondary protons. 

We consider an expanding universe homogeneously filled by the sources of accelerated 
primary UHE nuclei Aq with a generation rate per unit of comoving volume Qa (T, z) 
given by 

where 7 9 > 2 is the generation index, is the nucleon mass, and Cq is the source 
emissivity, i.e. the energy generated per unit of comoving volume and per unit time at 
z = 0. In Eq. ([55]) r min ~ 1 is assumed. In all calculations in this paper we assume 
also the existence of a maximum energy of acceleration E^ x = Zq x 10 21 eV (or = 
(Z q /Aq) x 10 12 ) with the condition Qa (T 9 ) = at T g > T^ x . 

Propagating in the space, a nucleus Aq experiences evolution, producing the secondary 
nuclei with different A and we calculate now their generation rate Qa(T,z) per unit of 
comoving volume. 

For calculation of the generation rates the evolution-trajectory formalism, developed 
in section 1231 is used. The trajectories A(z) = A(T, A, zq, z) and T(z) = Q(T,A,Zq,z) in 
the backward time allows us to calculate the generation parameters z g and r 9 , when A(z) 
reaches Aq. 

A{z) can be calculated as continuous quantity, but we use the method of ^4-jump 
trajectories, as it is described in section 12.21 Namely, the trajectories are calculated 
assuming in Eqs. (jSJ) for evolution of V A = const, until A reaches ^4 + 1 (or A + 2). 
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Transition A — > (A + 1) is assumed to occur instantaneously. Accordingly, we use in 
this case the jump behaviour of /3jL. In fact both methods, continuous A and jump A 
approximations, give the same results. 

Since a recoil momentum in the processes of transition, e.g. (A + 1) — >• A + iV, is 
negligibly small, one has approximate equality of Lorentz-factors of all three particles 

Ta+i ~ Ta ~ IV (56) 

The generation rate of secondary nuclei A can be found from conservation of number 
of particles. Consider an allowed trajectory which connects (Ta,za) state at generation 
A with (Tg, z g ) at generation of primary nucleus Aq. During time dt g at the epoch z g we 
have dNA primary nuclei generated with (T g , T g + dT g ): 

dN Ao = QA (T g ,z g )dT g dt g . (57) 

All these nuclei appear as dNA ^4-nuclei with (Ta, Ta + oITa) during time dtA- 
dNA = dNA ■ Thus, we have 

Qa(Ta, z)dT A dt = Q Ao (T g , z g )dT g dtg (58) 

where the states (Ta,z) and (T g ,z g ) are connected by the allowed trajectory. From Eq. 
([58]) and the time dilatation relation, 



dt g /dt = (l + z)/(l + z g ), (59) 
we obtain the generation rate of nuclei A as: 

QA(r A ,z) = Q Ao (r g ,z g )±±±§°. ( 60 ) 

i + z g ai a 

with dTg/dTA given in appendix [Bj 

As was emphasized above, the chain of equalities (j58[) is valid along the allowed tra- 
jectory, which breaks at some low A and z, where ta becomes larger than rp. If at z = z c 
7 + (A + 1) — > A + N is forbidden for all T < T c (z c ), then Qa(^a, za) = at za < z c and 
Fa < T c (za) (see also subsection 12.31 ). Starting from the critical redshift z c and below it 
^4-nuclei are not produced and (A + l)-nuclei propagate remaining undestroyed. 

Note, that for derivation we do not need stationary solution, and due to condition of 
homogeneity we do not need to discuss the spatial displacement of a particle. 

Coming finally to the generation rate of nucleons produced in (A + 1) — > j4 + iV 
transition, Qp +1 (T p , z), one may notice that this generation rate is the same that one 
for A nuclei, since the Lorentz factor of both particles are equal and they are produce 
simultaneously. 

Q* +1 {T a ,za) = QA (T g ,z g )\±^§!L, (61) 

i -+■ z g ai a 

where za and z g is the redshifts of A and Aq production, and Ta is the Lorentz factor of 
nucleus A and nucleon. We do not distinguish nucleon and proton because UHE neutron 
decays fast, on the time scale of this problem, to proton. 
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B Ratio of energy intervals at epochs of produc- 
tion and observation 



In this appendix we derive the ratio between Lorentz-factor intervals dT g /dT for a nu- 
cleus with variable atomic mass A(z) propagating along the evolution trajectory T g (z) = 
Q(A, T, zq,z), where the first three indices describe the initial conditions and z is the run- 
ning redshift. The meaning of the generation index g here is more general than the epoch 
of A{z) evolution to Aq, though it includes this case, too. T g {z) means here the Lorentz 
factor that a nucleus has at the running epoch of evolution z including z g . 

We use here the Lorentz-factor loss for a nucleus A due to pair production as b^ ah (T, z) = 
—dT/dt, and express it through 6g(T) = —dT/dt for proton at z = 0. At epoch z we have 

Z\z) t 



^pair 



A(z) 



(1 + z)%[{1 + z)T], 



(62) 



We introduce k = {Z 2 /A 2 ), assuming it be constant along the evolution trajectory A(z), 
and thus using Z 2 /A = kA(z). 

The Lorentz factor of nucleus A along the trajectory can be presented in a general 
form as 



T g {t 1 ,t 2 ) = T+ f 2 df 
Jti 



d£\ /dT 



dt'J 



ad 



d,t> 



pair 



Changing the variable t to z and using dt 

dz' 



T g (z ,z) = r + 



-dz/[(l + z)H(z)], we obtain 
,1 + z' 



7 r g (z>) + k I dz'— A(z') 6g[(l + z')T g (z% 



Differentiating Eq. (|64p in respect to T and using 

d6g(rp _ a&g(r') dV 
dr ~ or dr' 

one finds for the ratio of Lorentz-factor intervals y{z) = dT g {z)/dT: 

'dbl(T>)\ 



y{z) = 1 + 



dz' / /\ , t i , {l + z') 2 A{z') 
,y{z) + kjjz y{z) 



1 + z 



dT> 



T<=(l+z<)T a (z>) 

Differentiating Eq. ()66f) in respect to z we obtain a differential equation for y(z) 



1 dy(z) 



1 



+ 



k 



(l + z) 2 A(z) ( db%(T>y 



y(z) dz 1 + z ' H ^n m (l + z)* + A V dT ' 
The solution of Eq. (|67p can be easily found as 



(63) 
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r'=(i+*)r 9 («) 



dTg(z) (1+Z) 



dT{z ) (1 + Zq) 



cxp 



k_ f z , (l + z') 2 A(z r ) 

H Jz z ^/n m (i + z') :i + n A \ dT' 



T>=(l+zi)T g {z>) 

(68) 

Equation (168p gives the ratio dTg/dT, where dr is the Lorentz-factor interval at the 
initial state z$ (with most important case zq = 0) , when a nucleus has fixed A and Lorentz 
factor r, and dT g is the interval on the evolution trajectory at redshift z when the atomic 
number is A(z) and the Lorentz factor is T g {z). 

In the applications we often need dT g /dT ratio for a fixed A. The most important case 
is the primary nuclei Aq. The ratio for fixed A follows trivially from Eq. f|68[) substituting 
there kA(z') by constant value Z 2 /A. For the proton case this value is 1. 



O 1 



C Solution to the kinetic equation 

In this appendix we derive the solution of the kinetic equation for secondary nuclei from 
which the solutions for primary nuclei and secondary protons can be easily obtained. 

We consider the secondary nuclei A being produced homogeneously in the space with 
the rate Q A (T,t) and then propagate as a nucleus species with the fixed (unchanged) A 
until it is photo-disintegrated. The kinetic equation reads 



dn A (T A ,t) d n A (T,t) , , 

—di ^ A (r,t)n A (r,t)] + -^ = Q A (T,t), 



where b A (T,t) is the rate of Lorentz factor loss given by 

Mr,*) -- 



f = r^) + -^ pair (i>) 



(69) 



(70) 



with H(z) and 6^ air being the Hubble parameter at redshift z and pair-production loss 
for proton, respectively. The time of photo-disintegration is given by t^ 1 = Aj3^ is or by 
r^ 1 = dA/dt, see Eq. ©. 

The characteristic equation for the kinetic equation ([69]) reads 

dT/dt = -b A {T,t). (71) 

With T(t) taken on the characteristic in Eq. (j6"9|) . the term b A (T,t)dn A (T,t)/dT disap- 
pears and the kinetic equation (169|) takes the form 



dn A (T A ,t) 
dt 



+ n A (T,t) 



ffl&ir(r,i) db ad (T,t) 



or 



dr 



+ r A \T 1 t) 



Q A (T A ,t) (72) 



Using the notation 



Pi(z) = db ad (z)/dT = H(z) 

P 2 (z) ^ db^{T 1 z)/dT = ^{l + zfW)/dT') T , ={l+z)T , 



(73) 



and taking into account that with V on the characteristic, time t becomes the only variable, 
one obtains the solution of Eq. (|72j) as 



n A (t) 



dt'Q A (t')exp 



fdt"(-p x {t")-p 2 [t")+Ti\t": 



Changing the integration variable t to z and using 

dz 1 

dt = 



dz 



(1 + z)H(z) " H (l + z)^£l m {l + zY + fi A ' 
the solution can be written as 



(74) 



(75) 



exp 



dz" 



n A (T,z = 0) 
(1 + z")H(z" 



exp 


[L 







dz" - 



■ Q A [T'(F,z')} 
(1 + z')H(z') X 

W) 

1 + z")H(z") 





" r*o dt" ' 


exp 


[ U T A (T,t")\ 



(76) 



on 



We keep the last integration over t" to make clear the physical meaning of this integral 
as a suppression factor for the survival time of nucleus A. The upper limit to is the age 
of the Universe. 

Putting Pi(z") and P2(z") from Eq. (|73p into Eq. (|76p . we find that the product of 
the first two exponents gives the ratio of energy intervals calculated in appendix [Bl 



(1 + z')exp 



Z 2 1 



dz"- 



(1 + z 



ll\2 



a h J y/n m (i + z"f + n A \ dr 

Finally, we have 

n A (T,z = 0) 

where 



r"=(l+z")T a (z") 



dTf(z'^ 



dT A 



dz , Qa[t'(t,z')] dv e _ v{T ,^ 



V (T',z' 



to 



dt" 



1 + z')H(z') dT 
1 



dz"- 



1 



(77) 
(78) 

(79) 



{1 + z")H{z")t a {T",z") 
oo as upper limit, since it is regulated by the 



T A (T",t") 

In fact, in Eq. (J78"j) one can put z ma _ x 
factor exp[— r)(T', z')]. 

The equations (fTH|) and (jT9j) are valid for primary nuclei (Aq,Zq) and secondary pro- 
tons; in the latter case one should put t a — > oo and r\ — > 0. 



D Comparison of the secondary nuclei and sec- 
ondary protons fluxes 

We perform here the analytic comparison of the fluxes of secondary protons and secondary 
nuclei at the same Lorentz factor T in the most general form valid for kinetic equations 
and demonstrate that fluxes of the secondary protons are always higher. 

The secondary nucleus A and ^4-associating secondary nucleon N, born in the same 
decay (A + 1) — > A + N, are twin brothers: they are born at the same redshift z with the 
same Lorentz factor F and with the same rate of production Qa(T, z) = Qp(T, z) (we do 
not distinguish between neutron and proton because of the fast decay of the former). But 
one of the brothers, A- nucleus, lives short time ta and its flux is lower. We estimate here 
this effect analytically. 

The flux of the ^4-associating secondary protons is discussed in sections 13.21 and [3] and 
in the most general form is given as 



,max 



A (r m _ / - dz' n A+1 (T',z') ( dV\ 

where zo = is the initial condition, the generation term Qp(T,z) = Qa(F,z) is written 
in the most general form, and T'(z') = G p (T,0, z') is the proton trajectory, on which 
dT 1 /dT ratio is given by 



dT'(z') 
dT 



= (1 + z') exp 

p 



dz"(i + z ") 2 f db p (ty 

H(z") [ df 



f=(i+z")r"_ 



(81) 



oo 



Flux of the ^4-nuclei is calculated as 



with 



»u(r,o) 



dT'(z') 

dr 



max 
-A 



dz' 



n A+1 (T',z') 



(1 + z>)H(z')t a+1 (T>,z>) 



_____ 

dr 



-V(T'.Z>) 



(82) 



(1 + z')exp 



•-' dz"(i + z") 2 (<n%(fy 



H(z") 



f=(l+z")T" 



(83) 



where we assumed < Z/A >= 1/2, 

One may observe similarity of expressions for secondary protons and nuclei. It is easy 
to deduce inequality n p (T,0) > n A (T,0) caused by two main reasons. The first one is 
given by 2;™ ax being considerably larger than z^^, since the latter is limited by the decay 



lifetime T4 (see factor e _r? in Eq. I82[) , The second reason is connected with larger energy 
loss of nuclei to pair production /3a(T) = (A/A)f3 p (F). It results in T' A > T' at the same z' 
in the integrand of Eq. (|82p . and thus to suppression of this integral. The larger dV /dT 
for nuclei compensates partly the discussed flux inequality. Since this inequality is valid 
for every A, the total flux of the secondary protons also exceeds that of nuclei: 



£ A <(r,o) >s A n A (r,o). 



(84) 



At equal energies the nuclei flux should be shifted by factor A to higher energies and in 
practical cases it exceeds the secondary proton flux. 
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